#!/usr/bin/env python
"""
analysis function
"""

__author__ = 'Oliver Lindemann'
__version__ = '0.4'

import numpy as np
import math
import data
import scipy.stats.stats as stats
import sys

def log_regression(x , y):
    """ log regression between x & y
        log(y) = exp*log(x) + log(a)
            or 
        y = a * x^exp 
        returns exp, a, correlation (r2)
    """
    exp, fact, r2 = stats.linregress(np.log(x), np.log(y))
    return exp, np.exp(fact), r2

def SNARC_slopes(target, button, RT):
    """ calculates median RTs and SNARC slopes
        returns result-dict:
           beta, interc, r2, median_left, median_right, delta, 
            n_left, n_right, target_numbers, 
        numpy vectors:
            target: target numbers, RT: reaction times
            button: 1 left , 2 right , 0 no button
    """
    print "calculating SNARC slopes"
    target = np.array(target).squeeze()
    button = np.array(button).squeeze()
    RT = np.array(RT).squeeze()
    if target.shape != button.shape or button.shape != RT.shape:
        print "ERROR: unequal vector length"
        sys.exit()

    rtn = dict(beta=None, interc=None, r2=None, delta=None)
    rtn['target_numbers'] = data.get_values(target)
    for btn in [1, 2]:
        idx = np.nonzero(button == btn)[0]
        T = target[idx].copy()
        R = RT[idx].copy()
        M = np.array([])
        N = np.array([])
        for n in rtn['target_numbers']:
            idx = np.nonzero(T == n)[0]
            M = np.append(M, stats.nanmedian(R[idx]))
            N = np.append(N, len(R[idx]))
        if btn == 1:
            rtn['median_left'] = M
            rtn['n_left'] = N
        else:
            rtn['median_right'] = M
            rtn['n_right'] = N

    if len(rtn['median_left']) == len(rtn['median_right']):
        rtn['delta'] = rtn['median_right'] - rtn['median_left']
        rtn['beta'], rtn['interc'], rtn['r2'], prob, stderr = stats.linregress(rtn['target_numbers'], rtn['delta'])

    return rtn

def quintile(data_vector):
    """ returns vector of the same length 
        containing beloning quintile bins (0-4 & nan)
        
        function usefull for RT bin analysis
        example:
            q = quintile(rt_data)
            for x in range(0,5):
                idx = np.nonzero(q==x)
                print np.median(val[idx[0]])

    """
    val = data.remove_nan(data_vector.squeeze())
    val.sort()
    quint_size = (len(val)) / 5
    quint_threshold = []
    for x in range(1, 5):
        quint_threshold.append(val[x * quint_size])
    quint_threshold.append(val[-1] + 1) #larger than largest eleemnts

    #write quintiles 
    quint = np.array([])
    for x in  data_vector.squeeze():
        idx = np.nonzero(quint_threshold > x) #which threshold is lower 
        if len(idx[0]) > 0:
            quint = np.append(quint, idx[0][0])
        else:
            quint = np.append(quint, np.nan)

    return quint


def cramers_phi(chi2, N, length_of_minimum_dimension=2):
    """effect size for Chi squares"""
    df = length_of_minimum_dimension - 1
    ef = np.sqrt(chi2 / (N * df))
    print "chi2=", chi2, ", df=", df, ", N=", N , "::::: PHIc=", ef
    return ef


def tvalue_crit(df, alpha):
    """ returns critical t value (two-sided testing"""
    prob = 1- ((alpha)/float(2))
    return stats.distributions.t.ppf(prob, df)

    
def  repeated_measures_confidence_interval(MS, df, n, alpha=0.05):
    """returns the confidence interval (-+CI) for within-subject designs 
    as recommended by Loftus & Masson ('94)
    
    Note: A confidence interval of mean differences is CI*sqrt(2); not(CI*2) 
    
    function CI = rmConfidenceInterval( MS, df, n, alpha)
    <MS>: mean square error
    <df>: df of MS
    <n> : number of subjects
    
    returns: -+CI 
    example article: repeated_measures_confidence_interval(35.32,27, 10, 0.05)
    
    Oliver Lindemann"""

    prob = 1- ((alpha)/2)
    crit_t = stats.distributions.t.ppf(prob, df)

    CI = math.sqrt( MS / float(n) ) * crit_t

    print 'CI = {0} (MS={1}, df={2} per group n={3}, two-tailed p={4}, crit t={5})'\
            .format(CI, MS, df, n, prob, crit_t)
    return CI



def visual_angle(viewing_distance, stimulus_size):
    """ returnes the visual angle"""
    size = stimulus_size / 2.0
    va = 2.0* math.atan(size/viewing_distance)
    print "Visual angle (radians) = " , va
    print "Visual angle (degrees) = " , va*180/math.pi
    return va*180/math.pi
